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SUMMARY 


Using global interpolation functions (ClIFs). boundary element solutions are 
obtained for two-dimensional laminar flows. Two schemes are proposed for 
handling the convective terms. The first treats convection as a forcing func- 
tion. and converts the flow equations to pseudo-Poisson equations. In the 
second scheme, some convective effect is incorporated into the fundamental 
solution used in constructing the pertinent integral equations. The lid-driven 
cavity flow is selected as the benchmark problem. 


INTRODUCTION 

The boundary element method (BEM) has traditionally been applied to prob- 
lems governed by linear differential equations. At the core of the basic BEM 
computational process is the fundamental solution (also referred to as the 
free-space Green’s function) defined as the impulse response of the governing 
equation to a unit action. This fundamental solution is either too difficult 
or impossible to derive for practical nonlinear problems. Recently, with the 
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introduction of the so-called Dual Reciprocity techniques (see e.g, f Nardini 
Sc Brebbia [1982]: Brebbia it at.. [1991]: Partridge it ai. [1992]; Cheng et 
a/.. [199:}]: Lafe [1993]: Lafe Sc Cheng [1994]). the method is being proposed 
for certain classes of nonlinear problems. 

Using the Dual Reciprocity approach, a given problem is typically decom- 
posed into two parts - the linear and nonlinear portions. The solution to the 
linear portion is represented by a boundary integral whose kernel consists 
of the fundamental solution to the linear governing equation. The nonlin- 
ear part is represented by either 1) local bases functions (Brebbia et ai, 
[1991]): or 2) global interpolation functions (GIFS) (Lafe [1993]). In either 
case, the boundary integral expressions and interpolation functions contain 
coefficients whose values are to be determined by enforcing the boundary 
conditions. When the "direct BEM“ approach is followed the unknown coef- 
ficients are in essence the unknown physical variables (velocity components, 
pressure, temperature) of the problem. On the other hand, using the “indi- 
rect BEM" approach, the unknown are the weights/strengths of the boundary 
sources/dipoles and the local/global interpolating functions. The computa- 
tional intensity of the indirect approach is much less than for the direct. 

In this paper, we develop a GIF-based indirect BEM code for two- di- 
mensional steady-state incompressible Xavier-Stokes equation. Test results 
are shown for the lid-driven cavity problem. 


GOVERNING EQUATIONS 


The governing equations are: 


da dv 
0.r dy 
du du 

d.r dy 


dv dv 

u aI- + % = 


i dp p / d 2 u d 2 u\ 

pd.r p U-t' 2 dy 2 ) 

i dp p ( d 2 v d 2 v\ 

p dy p V d. v 2 dy 2 ) 


( 1 ) 

( 2 ) 

( 3 ) 


where (t/.r) are the velocity components in the r and y directions respec- 
tively. p is the pressure, p is the density, and // is the viscosity. Let 

A' = jt/L 
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V = y/L 
r = ,,/r 

V = r/r 

P = p/(pr 2 ) 

With these the governing e(|iiations become: 


or or op i fo 2 r o 2 r\ 

OX + OY ~ OX + P, \0X 2 + OY 2 ) 

OY 0\' OP I f0 2 \ ' 0 2 V\ 

' OX + OY ~ ~ OY + It VrAV 2 + OY* ) 

where (lie Reynold s Number R, — pT'L/p. 


( 4 ) 


(5) 


(G) 


BOUNDARY INTEGRAL EQUATIONS 


In order to convert the above into boundary integral equations two ap- 
proaches have been followed. In the first approach, the entire system of equa- 
tions is converted into an elliptic system, with the convective term wholly 
embedded in the right-hand-side forcing function. There is concern about 
the suitability of the elliptic system to adequately represent the convective 
forces at moderate to high Reynold’s number regimes. The second approach 
rectifies this through a more direct perturbation-based analysis which is more 
suited to capturing convective effects as the Reynold s number increases. 

Approach I 


The above equations are converted into an elliptic system: 


v 2 r = c, 

(U 

= f 2 

(8) 

v 2 p = /••, 

(9) 


where 


„ „ (,.or t .or op' 

F ' - ' yv +l W + 7Jx, 


235 


ORIGINAL PAGE IS 
OF POOR QUALITY 



i 


F 2 

F* 


Dl ,.c)V t .dV 

*- {l ax + x W 


A 


di’dv dr m 


\dxcjy or ox 



(10) 


Boundary Integral Equations 

The associated indirect boundary integral equations are: 


L'(X.Y) = 

/ ir,(A". Y')g(X'.Y’: X.Y) (IT + £ 3 lk Y) 

k 

(11) 

V(X.Y) = 

f r u' 2 (X'.Y’)g(X'.Y': X.Y) (IT + £ ^k(X. Y) 

(12) 

P(X.Y) = 

[ w 3 (X'. Y')g(X'. Y': X. Y)(IT + Y A', Y) 

k 

(13) 

which 

Q — In [(A' - A") 2 + (V - V') 2 ] 



T:i lk M,(X.Y) = F,(X.Y) 

(14) 


k 


X !i V k (X.Y) = M k (X.Y) 

The functions M k { X. Y) are t he interpolation functions used in represent- 
ing t he convective terms. If we choose 


Mk = co#(m k X)coM n k ) ) 
it is easily shown (Lafe [1993]) that 

_ co#(m k X)cox\n k Y) 

* k , ) i , 

(nf. + mj : ) 


(15) 

(16) 


Approach II 

Our aim here is to have a better incorporation of the convective effects in 
the driving differential operator. Let 

r = r„ + r, 
v = h> + r, 

P = Po + Pl 
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where (f 0 . 1 o- Po) denote solutions to a convection-tree flow Held, while (U\, Vi, P\) 
represent the convective effects. Hence 


v 2 r„ 

rn 0 


V 2 P 0 


0 

0 

.J 01 o ( P o 

“ ox or 


()l 0 0\ 0 

W ox 


(17) 

(18) 

(19) 


while 


A.^r l -o l a\ ) .Vo.fXA\.p l ) = j\u\>.Vo.Po) + ih(Px.Vi) ( 20 ) 

•* r 

-^-V 2 lj -o-iH'o. I o./'|. I ,. P x ) = fAr v .V 0 .P 0 ) + h 2 (l\.\\) (21) 
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0) 
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The external boundary conditions are imposed on ( ( o- 1 o- Po)- Therefore the 
variables (l \ . \ \ . P\ ) are allowed to enjoy homogeneous boundary conditions. 
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Boundary Integral Equations 

The velocity components (l o. V 0 ), which are governed by Laplace's equa- 
tions. can be represented by ‘pure boundary integral equation using the 
fundamental solution for potential flow. However the pressure term, because 
of the non-zero forcing function, will include global interpolation functions. 

Therefore, the indirect boundary integral equations for the convection- 
free variables are: 

(MX. Y) = J w 0 i(X'. Y")g(X'. V: X. Y) (IT (23) 

lo(.V.V) = J r W 02 (X'.Y')g(X'.Y':X.Y)cir (24) 

Po(X.Y) = f w Q3 (X\Y') g (X\Y':XX ) f ir + Y. j o^^(X.Y) (25) 

Jr k 

where ( «' 01 . ic Q2 . u' 03 ) are fictitious sources while 4'*. are t lie GIFS. 

The convective effects (l ]. Vj./ 5 ,) are represented by GIFS. For flows in 
simple geometries it is possible to select GIFS which automatically satisfy 
the required homogeneous boundary conditions. 

Hence 


1\{X.Y) = Y.Ju*k(X.Y) 

k 

(26) 

Vi(-Y.V) = 5>*¥*(A\r) 

k 

(27) 

Pi(X.Y) = Xy^G-Y.r) 

(28) 


k 


in which (.in-, 3 2 k, h k) the pertinent coefficients for the GIFS. These 

parameters are calculated by enforcing 

1. Equations (20-22) at selected collocation points within the flow region 
and 

2. Homogeneous conditions at selected boundary points. 

The chief advantage of the first approach is the simplicity and size of 
the global coefficient matrices, derived from the Laplace equation solver. 
This translates into a compact, last, and highly efficient numerical imple- 
mentation. The drawback is its iterative character since the forcing function 


238 


ORIGINAL PAGE IS 
OF POOR QUALITY 



depends on the solution being sought. Convergence is difficult to achieve for 
large Reynold's number because the governing elliptic system becomes singu- 
lar. and cannot adequately represent the underlying physics of the problem. 
The second approach, which in essence' separates the convection-free flow 
from the main flow, allows for a more direct representation of the asymptotic 
limits of the Reynold's number. Furthermore, by products of higher-order 
terms (i.e., setting hi = h 2 = hi = 0. the solutions can be obtained without 
iteration. However, the coefficient matrix is larger and the approach involves 
a greater level of computational intensity. 


NUMERICAL IMPLEMENTATION 


Approach I 

We subdivide the boundary into n* elements. Let A*(x) (k = 1,2, 
represent the bases functions describing the distribution of w on T. In the 
examples being reported in this paper, constant bases functions are being 
used for the fictitious strengths tr, on the boundary. By selecting each of the 
/?& boundary points as successive origins of integration, the pertinent integral 
equations can be assembled into the system: 


= !> ' * - 1 

A = 1 

where 

f rk .V; (x')</(x'.x, ) c/x' X, € r* 

fr k V ( x' )0g/dn ( x'. x, ) fix' x, e Tq 

_ f <t>(x, ) - 3;' Vij X, G T.J, 

' “ \ d*/dn{Xi) - E."=, .ijMij/On X, € r Q 

where 4> = (f'.WP). Therefore, we have in, equations to determine uy 
(k = 1.2. • • • ni,). Symbolically equation (29) can be written in the alternative 
form: 

.4VV = B (32) 

which can be inverted to give: 



(29) 

(30) 

(31) 


W = A~ l B 


(33) 
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The whole process boils down to the iterative solution of equations (14) 
and (33), with repealed updating of F using (10). The iterative steps are: 

1. Start with a trial F (/'.(.. F, values for i = 1. 2. • • • n,y). 

2. Obtain 3 from equation (14). 

3. Obtain W using equation (33). 

4. 1. se discretized forms ol the appropriate integral equations to compute 
4>. V4> at all n d points. This provides a better estimate for F. 

5. Go back to Step 2 if convergence condition is still unsatisfied. 

Note that the matrix inversions in equations (11) and (33) need only be per- 
formed once, for fixed boundary problems. The vectors W and 3 are the 
quantities whose values change during the iterative process. Once conver- 
gence is reached, the discretized integral equations can be used routinely to 
obtain $ = ((.!'. P) or the gradient at any point (x) of interest. 

Approach II 

The numerical implementation for the convection- free quantities (f ’ 0 ,Vo,P 0 ) 
is similar to the one followed in Approach I. with t he coefficients for the GIFS 
set to zero for the velocities. No iteration is required. 

The convective- flow quant it ies are calculated t h rough the coefficients ( 3 \ i< , 82 k- • 3j<.) 
whose values are obtained by solving the following coefficient matrices: 

”T "T " T 

XI Jik^iik + XI + X! = F\< i = 1.2. ••• n d (34) 

A-l k-= 1 k= 1 

nj- tij nj 

X hk-B\,k + XI -huB-ak + XI hk B-^k- — F ,2 i = 1.2, • ■ • n d (35) 

*•= 1 A=1 k= 1 

n T 71 t 7 

XZ + X ^3 i-C-iik + X likC'Mk = Fn i = 1.2. ■ • • n d (36) 

A = 1 k = 1 k=l 

"T 

XI hkBtjk =0 j = 1 . 2. - • • m, (37) 

A = 1 
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Y2 Ou-Ei.it,- =0 j = 1.2. ■••»»» 

k- 1 

nr 

(38) 

YL ^lWi./A =0 j = 1.2. ■ • • /?;, 

A-l 

(39) 


where (if the higher order terms are neglected) 
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^tk 

B\ ik 
B’lik 
B*k 


Cuk 


( 2 ik 
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0 A < J) 0) 
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( x, ) 
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if x, efq, 


In the above <I> = ((*. V. P): Q = (01 ’ /On. 0\ /On . UP/ On): n j = 3(n<, + rid), 
and X, = (AW) for 2D flows. 
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TEST RESULTS 


We examined the lid-driven cavity flow problem depicted in Fig. 1. A unit 
horizontal velocity is imposed on the lid (at 1 = 1). while the no-slip bound- 
ary condition l = \ = 0 is imposed on all solid walls. The boundary 

condition for the pressure on all walls is (Fletcher [1991]): 

OP (0l_ _ 0V\ 

dn ~ R f Os VdV ~ OX ) 

A typical convergence profile, using Approach I. is shown ( R f = 15) in Fig. 
2. The horizontal velocity at the vertical center-line is shown in Fig. 3. 
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CONCLUSIONS 


A boundary element code, based on the use of global interpolation func- 
tions. for solving the Xavier Stokes equations have been proposed in this 
paper. The avoidance of any domain integration shows the enormous power 
of the technique. As long as the underlying physics of the problem is ade- 
quately represented in the fundamental solutions used as the kernel of the in- 
tegral equations, accurate simulations can be carried out for moderate to high 
Reynold s number flows. Only trigonometric bases have been used to repre- 
sent the nonlinear convective terms. Investigations are currently underway 
for employing other bases including those derived from orthogonal functions 
such Chebychev polynomials, wavelets, and cellular automata transforms. 
Three-dimensional OIF-based BLM code for internal flows are also being 
developed. 
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Figure 1: Lid-driven Cavity Problem 
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Figure 2: Typical Convergence Profile 
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